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SUMMARY 

We develop a multilevel scheme for solving the balanced long transportation problem, that is, 
given a set {c^j} of shipping costs from a set of M supply nodes 5* to a set of N demand nodes Dj, 
we seek to find a set of flows, \xkj}, that minimizes the total cost XkjCkj ■ We require that 
the problem be balanced, that is, the total demand must equal the total supply. Solution 
techniques for this problem are well known from optimization and linear programming. We examine 
this problem, however, in order to develop principles that can then be applied to more intractible 
problems of optimization. 

We develop a multigrid scheme for solving the problem, defining the grids, relaxation, and 
intergrid operators. Numerical experimentation shows that this line of research may prove fruitful. 
Further research directions are suggested. 

INTRODUCTION 


The transportation problem is the simplest of network flow problems. It is posed on a bipartite 
graph, consisting of a set of M supply nodes, a set of N demand nodes, and a set of arcs connecting 
them. Each supply node S% has a fixed amount s, of a commodity which it can provide. Each 
demand node Dj has a fixed requirement dj for that commodity, and for each arc (i,j) connecting 
supply node S { to demand node Dj there is an associated cost per unit flow Cij. When the total 
supply equals the total demand the problem is balanced. When M « N, the problem is referred to 
as a long transportation problem. Denoting the flow on arc (i, j) by Xij, the transportation problem 
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can be expressed 


M N N M 

Minimize subject to: J^Xjj = s ». YL x 'i = x *j — 0. 

*=i j=i ~j=i ^ »=i 

Let 6 denote an (M + iV)- vector whose first Kf entries are the available supplies at nodes Si 

through Sm, and whose last N entries are (negatives of) the required demands dj at demand nodes 

D\ through D^. Let K be the number of arcs in the problem. Throughout this work we shall 

assume that every supply node is connected to every demand node, so that K = MN. Let the 

if-vector x be composed of the flow on the arcs from the M supply-nodes to the N demand nodes 

in some order, and the K- vector c be the cost of shipping on those arcs in the same order. We j 

denote by A the incidence matrix of the graph, so that A has as many rows as there are nodes in 

the problem, M + N, and as many columns as there are arcs (MN). Each column of A is 

associated with one arc of the problem, and they are arranged in an order that matches the order of 

the vectors c and x. Each column has exactly two non-zero entries: a +1 in the row corresponding 5 

to the tail (supply) node Si of the arc, and a —1 in the row corresponding to the head (demand) I 

node Dj. Each row of A is associated with one of the constraints of the problem [ 1 ]. Then the » 

problem may be written in matrix notation as * 

Minimize: c T x - 

Subject to: Ax = b, 

x > 0. i 

~.m, r . • _ . _ ’ ' i 1 

A simple example is presented in Figure 1. In the example, there are four supply nodes, having g 

12, 15, 10, and 7 units of the commodity to deliver. There are three demand nodes, requiring 13, 

20, and 11 units of the commodity. We seek to find a flow vector 1 

( \T I 

^ Xu X\2 Zi3 X21 X22 ^23 %31 %32 ^33 ^41 ^42 ^43 J 

given that tTie vector of costs, written in corresponding order, is 

(215643174234 ) T . I 

The algebraic description of this problem is to find x such that c r x is minimized, subject to the 7 

system of constraints 

ft ;) ■ 
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Figure 1: A simple example of a transportation problem 

Very little work has been done on multigrid methods for discrete optimization problems. 
Significant studies to date are [2], [3], [4], and [5]. The traditional optimization algorithm which 
most closely resembles a multilevel algorithm is aggregation/disaggregation [6], [7], and [8], in 
which nodes are aggregated in a logical way in order to reduce the size of the problem, and the 
solution to the smaller problem is disaggregated to provide an initial estimate for the solution to 
the original problem. The most successful work to date, and the work that inspired this study, is 
that of Kaminsky [4] . 


COST-SPACE 


In [4] it is required that the demand nodes occupy a physical location in space, and that a 
relationship exist between transportation costs and distances. This is done so that the coarsening 
step may be performed by aggregating together demand nodes that are physically near one another. 
For this to make sense, it is necessary that shipment to each of the aggregated demand nodes 
involve a similar cost, which naturally occurs if the shipping cost is a function of distance. For 
many applications this makes perfect sense; the cost of shipping a commodity is often directly 
linked to the distance the commodity must be shipped. This restriction is overly limiting for other 
types of problems, however. For example, the manpower assignment problem, in which a specified 
number of jobs must be assigned a given set of workers, can be formulated as a transportation 
problem. There is no distance involved in such a problem, and cost of assignment is related to other 
factors, such as the cost of training an individual for a specific task. 

In order to address problems that have no geometrical dependence of cost on distance, we 
employ a change of coordinate systems from physical space to a space we describe as cost space. 

For the M x N problem, cost space is the M-dimensional space in which each of the coordinate 
axes is the cost of shipping from one of M supply nodes. Each of the N demand nodes is placed in 
cost space at the point whose coordinates are the unit costs of shipping from the supply nodes to it. 
For example, the three demand nodes in Figure 1 would be placed in a four-dimensional cost space, 
and would have the coordinates A = (2, 6, 1,2), A = (1,4, 7, 3), and D$ = (5, 3, 4, 4). This change 
of coordinate systems means shipping cost becomes the metric of the problem, so that two demand 
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nodes are “near” each other only if the shipping costs are similar, and the aggregation of 
neighboring demand nodes automatically ensures the similarity of their costs. 

Posed in cost space, the dimensionality of the problem equals the number of supply nodes. In 
traditional multigrid methods, one typically uses grids that are tensor products of one-dimensional 
grids, each having a cardinality of gridpoints that is a power of two. In the cost space approach this 
would lead to a very rapid growth in the size of the problem; for this reason the cost space 
approach can be applied only to problems with a relatively small number of supply nodes. This is 
one reason for restricting our attention to the long transportation problem. 

Reduced dimension cost space 

If at least one supply node is connected to all demand nodes (and in our work we assume this to 
be true of all supply nodes) then we can transform the M x N transportation to an equivalent 
(M — 1) x N problem, which we call the reduced dimension problem. Since we are dealing with the 
long problem, the transformed problem is somewhat simpler and less expensive to solve. The 
transformation is accomplished as follows. Suppose that supply node Si is connected to all demand 
nodes. Then for each demand node Dj, we subtract cy, the cost of shipping from supply node Si to 
Dj, from all of the shipping costs into demand node Dj. That is, we form an auxilliary cost vector 
Cij = Cij — cij. The result is that for supply node Si, all the demand nodes map to the origin in cost'* 
space. Effectively Si has been removed from the problem, leaving an (M -l)xiV problem to be 
solved. For example, if we use the cost of shipping from S 2 on the example in Figure 1, the 
transformed cost vector becomes 

c = ( -4 -3 2 0 0 0 -5 3 1 -4 -1 3 ) T . 

We can show that while the objective function value is different for the new problem, a solution for 
one is equivalent to a solution for the other. 

Theorem 1 Let the M x N balanced long transportation problem be represented by a bipartite 
graph G, and suppose that supply node Si is connected to all demand nodes. Let b b< the ( M + N ) 
length column vector whose first M entries are the supplies at the supply nodes and whose 
remaining N entries are the negatives of the demands at the demand nodes. Let A be the adjacency 
matrix of the graph G; that is, for each arc (i, j) we have A(i , (i — 1 )N + j) = 1 and 
A(M + j , (i-l)N + j) = -1. Let c be the (M + N) length vector whose k — (i — 1)N + j element is 
the cost Cij of shipping from node S, to node Dj along arc ( i,j ). Define c to be the vector whose k th 
entry is c* = dj — cij. Then x* is a solution to the problem 

Minimize 
Subject to 

if and only if it is a solution to the problem: 

Minimize: c T x 

Subject to: Ax = 6, 

x > 0 . 


Ax = b , 
x > 0 , 


64 


I! Rllll 11 I III IIII'II'1 1 1 !! II l Mill 1 mm t r » i ■miiniiri n i 1 1 


Proof: 


- T 
C 1 X 


M N M N 

= 53 E &kj x kj — E E ( Cfc J _ C i]) X kj 

fc=l j= 1 fc=l j= 1 

M N 

— 53 53 ( c kj x kj ~ CijXkj ) 

fc=l J=1 

MW MW 

= EE Ckj x kj ) 1 53 C h X kj 

Jfc=l j=l fc=l j= 1 

W M W 

= c T x - 53 Qj 53 x v = c T x - E c o d J» 

j=i fc=i i=i 


M 

since 53 x *j = in the balanced problem. But cydj does not depend on x, and therefore 
fc=1 

c T x achieves its extreme values precisely when c T x does. I 

This transformation of the costs to reduced-dimension space maps all costs of shipping from Si 
to the origin in cost space. As will be shown in the next section, our algorithm requires that the 
demand nodes be sorted once according to the cost of shipping. Since sorting is a fairly expensive 
operation, the savings generated by reducing the dimension of the problem are tangible. Once the 
transformation to reduced-dimension cost-space has been performed, the resulting problem may be 
solved with no further consideration of the transformation. Therefore, in the remainder of this work 
it is assumed that when an M x N problem is to be solved, it may be the reduced dimension 
version of a problem that was originally (M + 1) x N. 


A MULTIGRID APPROACH TO THE TRANSPORTATION PROBLEM 


Following traditional multigrid design approaches, we develop the necessary tools to devise a 
multigrid U-cycle, which we will combine with a nested iteration to create an FMG algorithm. In 
particular, it is necessary to devise restriction and prolongation methods, some form of local 
relaxation, and to weave them into an algorithm. 

Restriction 


To devise a restriction algorithm, it is first necessary to define a coarse grid. We use an approach 
in which each gridpoint on the coarse grid is a demand node for the coarse grid problem, and 
represents a pair of demand nodes on the fine grid. This is accomplished as follows. The demand 
nodes are first sorted by increasing cost of shipping from Si, and divided into two groups about the 
median of the sorted cost. This procedure results in two groups of demand nodes, one with a lower 
cost of shipping from Si, and one for which shipping from S\ is more expensive. Each of these 
groups are then sorted according to increasing cost of shipping from S 2 and divided into two groups 
about the median cost. This results in four groups, one for which shipping is expensive from both 
supply nodes, one group for which shipping is inexpensive from both supply nodes, one group for 
which shipping is expensive from S 2 and inexpensive from Si, and one group where shipping is 
expensive from Si and inexpensive from Si- 
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Figure 2: A simple example of a coarsening process for the transportation problem 


If there are more than two supply nodes in the problem, this process is continued. The four 
groups are each sorted by cost of shipping from supply node S3, then divided into smaller groups if 
necessary and sorted again, according to cost from St, etc. If the groups contain more than two 
nodes after the nodes have been sorted according to cost from all supply nodes, the sorting process 
begins again with cost from S\ on each of the groups. Eventually, the nodes will be sorted into 
pairs that have similar shipping costs from all supply nodes. Each of these pairs of demand nodes is 
then replaced with a single coarse grid demand node, the collection of which constitutes the first 
coarse grid. 

Further coarsening is accomplished by repeating the procedure described above on the coarse 
grids to produce still coarser grids. Figure 2 shows a simple example of the coarsening process. If 
the number of points on the original grid is a power of two, then in the limit a coarsest grid would 
consist of a single demand node. As in traditional multigrid methods, once the hierarchy of grids is 
established it is stored, so that the sorting process need never be repeated. 

Three quantities must be restricted when aggregating a pair of fine grid demand nodes into a 
coarse grid demand node: the demands, flows, and costs. Let D fff be the coarse grid node 
representing the fine grid nodes D’f and . It seems natural that the demands can be restricted 
simply by summing the demands of the two fine grid nodes to produce the demand at the coarse 
grid node, = Iff 1 [d),df{] = d* + dff . Similarly, the flow x 2 ^ from any supply node S k into the 
coarse demand node should be the sum of the flows from S k to each of the fine grid demand nodes 
that make up the coarse grid node, = x k j + x k i- 

Restricting the cost of shipment is more complicated, and no obvious “best” approach is 
apparent. However several methods can be considered. The simplest of these is to define the coarse 
cost to be the minimum of the fine costs, i.e., cjf^ = /^[c^, c£j] = min(c^ J , c%f). Other simple 
schemes are readily devised, such as using the maximum of the fine costs, or a weighted average of 
the fine grid costs. We use a weighted average of the fine grid costs. Again, there are several 
possible weightings, each having valid arguments for and against it. Three schemes were tested in 


some depth, equal weighting, flow weighting, and demand weighting: 


Equal weighting: 

^km 


Demand weighting: 

r 2h 

C'km 

II 

Flow weighting: 

h 

^ km 

= 


c kj + c k i 
2 ’ 

djC h kj + df c& 

d) + 4 ’ 

X kj C kj X kl C kl 

X h | 

fc.? ^ 


With flow weighting, provision must be made for the case where there is zero flow on both arcs. In 
such a case flow weighting can be replaced with either demand weighting or equal weighting. In 
general, we found that demand weighting most consistently gave the best results, and adopted it for 
our algorithm. 


Prolongation, or Interpolation 


Suppose that the problem has been solved on the coarse grid Q 2h . We seek a method of 
prolongation, that is, a way in which the coarse grid solution can be interpolated onto the fine grid. 
In the coarse grid solution there is some quantity of flow giving the flow from each supply node 
Sk to each coarse grid demand node d^. Each such demand node on the coarse grid, however, 
represents the aggregation of two demand nodes on the fine grid, d * and df . An interpolation of the 
coarse grid solution, therefore, can be constructed by treating the M flows . . . ,x\ J m , into 

the coarse grid demand node d^ 1 , as supplies. Interpolation, then, consists of solving for each coarse 
grid demand node, the M x 2 transportation problem with those supply values, the two demand 
nodes dh and df, and the shipping costs c k — 1 , 2 ,..., M. (Figure 4 shows schematically how 
the interpolation process appears.) 

Having defined the interpolation process as finding the solutions to many small transportation 
problems, we turn our attention to the mechanism for finding these solutions. A method for solving 
such Af x 2 problems is described below. The method is a special case of Vogel’s approximation 
method. 


Algorithm 1 Solving the M x 2 Balanced Transportation Problem 

1. For each supply nod( Sk- find the difference in cost of shipping 6k = | c^- — c£ ( | to the two fine 
grid demand nodes dh and d 1 } . 

2. Rank the M supply nodes in decreasing order of these cost differentials, so that 
6\ >62 > ... > S M ■ 

3. Repeat until all supply nodes are removed from the problem: 

(a) Denote the supply node eit the top of the ordered list as the “ current " supply node, and 
allocate flow to the demand node with the lower cost of shipme nt . that is, along the feast 
expensive arc, thus determining a “current" demand node. (In the event that more than 
one node has the largest differential cost, select from among them the node with the 
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<5 S = 1 


C 5 1 = 4 
C52 — 3 


Figure 3: Example problem illustrating the solution method for an M X 2 problem. 


smallest cost along one of its two arcs). Allocate flow along this arc until either the 
demand at the current demand node is satisfied or until the supply at the current supply 
node is exhausted. 

(b) If the supply at the current supply node is exhausted, remove that supply node from the 
problem. 

(c) If the demand at the current demand node is satisfied, remove that demand node from 
the problem, allocate the remaining supply from the current supply node to the remaining 
demand node, and remove the current supply node from the problem. 


Stop. 

As an example of this procedure, consider the five by two problem shown in Figure 3. The five 
supply nodes Si, S 2 , . • - ,S 5 have, respectively, 15, 12, 16, 18 and 14 units of the commodity to 
deliver. The demands of the two demand nodes Di and D 2 are 30 and 45. Let 6 = (4 8 4 6 1) T be 
the vector whose i th entry is the difference 6, between shipping cost from supply node S { to the two 
demand nodes (the costs themselves are given for each arc in the figure). Sorting from largest to 
smallest value of S i7 the supply nodes are ordered (S 2 , Si, Si, S3, S 3 ). Note that, while the 
differences for nodes Si and S 3 are the same, the cost cu along the arc from node Si to node D 2 is 
less expensive than either of the arcs incident from node 53 . Starting with node S 2 , then, as much 
flow as possible is sent along the least expensive arc. In this case, that is the arc to demand node 
D 2 . Since this demand exceeds the available supply from node S 2 , all of the flow from node S 2 goes 
along this arc. Similarly, node S 4 and then node Si send all of their supply to node D 2 . When node 
S 3 has sent 12 units of flow along its least expensive arc, the demand at node D\ is completely met. 
Thus node S 3 sends its remaining units to node D 2 , as does node S 3 . Although the arc from node 
S 5 to Di is less expensive, the demand at D ] has been met from supply nodes where the difference 
in arc costs is greater. 

We can show now that because of the special structure of the M x 2 problem, i.e., the fact that 
there are only two demand nodes, this procedure produces an optimal solution. 
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Theorem 2 Let x be the vector of flows assigned for the M x 2 problem using the algorithm given 
above. Then x is an optimal solution to the M x 2 problem. 

Proof: Suppose that x is not an optimal solution. Then there exists a flow x* ± x such that 
z* = c T x* is optimal. We will show that if x is determined by the algorithm given above and 
z = c T x, then z* > z, contradicting the assumption that x is not an optimal solution. Letting 
S k = | Cfc i - c k2 1 for each k = 1, 2, m . . . , M, assume the supply nodes have been ordered in 
decreasing order of 6 k so that h\> S 2 > ... > S M . Let i be the first supply node for which x* differs 
from x, and without loss of generality, assume that < c i2 . Let A = Xu — x* r We first observe 
five useful facts: 

1. Since the problem is balanced, total flow out of S k equals the supply, so that 
s k = Xfci + x k2 = x kl + x* k2 for every k, implying x fel - x* kl = x* k2 - x k2 . 

2. In particular, since A = x*i — x'h then —A = x\ 2 — X&. 

3. Since the problem is balanced, total flow into A equals demand, so that d x = £jii Xj\ and 
d\ — x l \ • Subtracting these two relations, and noting that x and x* do not differ for 

j = 1, 2, .... i — 1, we find that 0 = x n - x^ + Ejli+iObi - x* a ), implying that 

A = EjS l+ i(a5 1 -*ii). 

4. Using similar reasoning, we obtain —A = £*i i+1 (x* 2 — Xj 2 ). 

5. By construction, since c*i < c i2 , then Xn is as large as it can possibly be, so that if x and x* 
differ for node i then x« > x* lt implying that A > 0. 

Next, we observe that 

z* = z + z* — z 

= z + c T x* — c T x 

M M 

= z + 52 (xjiCji 4- x j 2 Cji) — 52 “f" x j 2 c j 2 ) ? 

J=i j-i 

where we have used the fact that x* and x do not differ for j < i. Separating the flows from supply 
node i, we can write 

M M 

z’ = Z + (xjj - Xn) Cil + (x- 2 - X i2 ) Ci2 + ( X }x “ X Jl) + 52 \ X h - X 3 2 ) C J 2 

j=i+l J=i4-1 

M 

= z A Cji + A Cj 2 + y ^ {^ji ( c ;i ^ 2 ) 

j=i+l 

where we have used the fact that x k \ — xj^ = x* k7 — x ki f or all k. Recalling that 6 k = |cfci — c^l, 
and that since Cn < C& then 6i — cn — c i2 , we observe that 

M 

z* > z + A Si + 52 ( x ji 

J=*+l 
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Figure 4: 

An illustration of interpolation by local optimization. The flows in the 3x2 coarse grid problem are 
the supplies for the two 3x2 local problems. The combination of the flows solving those two local 
problems makes up the interpolated solution to the 3x4 fine grid problem. 


Since the nodes are ordered in decreasing order of 6k , we know that — < — 6j for all j > i, and 
therefore 

M 

z* > z + A6i - Si Y ( x h ~ x n) ■ 

j=n+ 1 

Finally, recalling that A = Y!>jLi+i( x ji — Xji), we obtain 

z* > z + A Si — Si A. 


Therefore z* > z, contradicting the assumption that x is not an optimal solution. 


I 


Relaxation 


Suppose that we have solved the coarse grid problem, and have interpolated that solution by 
solving an M x 2 transportation problem for each coarse grid demand node in order to pass the 
local solution to the two fine grid demand nodes represented by each coarse grid node. The supplies 
for this local M x 2 problem are the coarse grid flows. Figure 4 displays a schematic showing how 
the interpolation process appears graphically. 

It is important to note that while each of the local M x 2 problems has been solved optimally, 
there is no reason to expect that the total set of fine grid flows thus assigned will be optimal. For 
this reason, it is essential that we devise some kind of “relaxation” scheme, whose task is to smooth 
or correct errors left by the interpolation scheme. 
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When two M x 2 subproblem solutions are viewed from a more global perspective, as a solution 
to an M x 4 problem (or as a portion of a solution to a still larger problem), this combination of 
locally optimized solutions may be flawed, in that too many arcs may have flow on them. This is 
because the minimum value of the objective function for a balanced transportation problem can 
always be obtained with a flow regime having flow on at most M + TV — 1 arcs. This simply reflects 
the fact from linear programming theory that an extreme point solution has flow on M 4- TV — 1 
arcs, if the solution is non-degenerate [9], and that an optimal solution can always be found at one 
of the extreme points (a degenerate solution is one in which distinct subsets of demand nodes are 
supplied by distinct subsets of supply nodes) . If the solution is degenerate, there will be fewer arcs 
with flow on them. For example, if TV > M, then in the extreme degenerate case each supply node 
provides flow to a disjoint subset of the demand nodes. This means that each demand node has 
exactly one arc with flow on to it, so that precisely TV arcs have flow. If TV < M, then the extreme 
degenerate case is when each supply node has exactly one arc with flow, giving M such arcs. 

When interpolating from Q? h to each coarse demand node generates two fine grid demand 
nodes and the optimal solution to the M x 2 subproblem has M -\- 1 arcs with flow, in the 
non-degenerate case. If the subproblem solution is degenerate, then M arcs have flow. If there are 
TV/2 demand nodes on Q? h , then after the interpolation the collection of subproblem solutions 
(viewed as the initial feasible solution to the Q h problem) will have flow on at least NM/2 and at 
most (TVM + TV)/2 arcs, depending on how many subproblems are degenerate. 

Thus, whenever NM/2 is greater than M + TV — 1, (which is true for any long transportation 
problem where M > 2 and TV > 3), the collection of local solutions has too many arcs with flow to 
be an extreme point solution to the fine grid problem, and is probably less than optimal. The local 
relaxation scheme developed here is designed to reduce the number of arcs with flow for the fine 
grid problem, which will generally have the effect of moving the global solution toward an optimal 
solution. 

The mechanism by which we do this is cycle removal. Since there are M + TV — 1 arcs in a 
spanning tree over M + TV nodes, and the addition of a single arc (or more) to a tree results in a 
graph with at least one cycle, then for most problems, the interpolation process will introduce 
cycles. We note that while this has been developed in the setting of the entire collection of local 
solutions, it is also true in a pairwise sense. That is, each of two M x 2 local solutions will have 
either M + 1 or M arcs with flow. Viewing the pair as a solution to an M x 4 problem, we observe 
that the combined solution will have at least 2 M arcs with flow. If TVf > 2 this equals or exceeds 
the M + 3 arcs with flow that would be present in an extreme point solution. 

To illustrate this, consider the possibilities when two 3x2 local solutions are combined into a 
3x4 solution, as shown in Figure 5. In a), two non-degenerate solutions are combined. Numbering 
the demand nodes of the combined problem clockwise from the upper left and the supply nodes 
from top to bottom, we observe that there are three cycles in the combined solution 
(Si, D 3 , S 2 ,D 1 , Si), (Si, D 3 , S 3 , D a , s 2 , Di,Si), and (S 2 , D 3 , S 3 ,D 4 , S 2 ). In b), a degenerate solution 
is combined with a non-degenerate solution, yielding a combined solution with one cycle. In c), two 
degenerate solutions are combined into a solution that has no cycles, while in d) , two degenerate 
solutions are combined into a solution that has one cycle. 

A reasonable candidate for a local relaxation process is to adjust the flow in the initial solution 
produced by the interpolation process, so that cycles are removed and the objective function is 
reduced. The effect of this procedure is to adjust the locally optimal flows which result from 
interpolation so that they are more nearly optimal in the global problem. 
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Figure 5: Combining two 3x2 solutions into a 3 x 4 solution. Four typical cases arc shown. 

Cycles are detected in the algorithm using a depth first search (DFS). The DFS proceeds as 
follows: 

Algorithm 2 Depth First Search 

1. Initialize all nodes with DFS number 0, to indicate they have not yet been visited. 

2. Start at any node. Assign this node a DFS number of 1, and define node 0 to be the 
predecessor of th is node. 

3. If any node adjacent to the current node has been previously visited , and has a DFS number 
lower than the predecessor of the current node , then the path from that node through the 
current node and back is a cycle. Stop DFS and call the cycle removal routine. 

f. If no adjacent nodes have lower DFS numbers, then look for any adjacent nodes which have 
not been visited. If there are any u n vis it eel adjacent nodes, identify the current node as the 
predecessor of the unvisited node, make the unvisited node the current node, and assign the 
current node a DFS number equal to the DFS number of its predecessor plus 1. 

5. If there are no unvisited nodes adjacent to the current node, make the predecessor of the 

curre nt node the current node. If the current node is node 0, stop. Otherwise, return to step 3. 

Once a cycle is detected, a cycle removal algorithm is used to adjust the flows. The technique is 
illustrated in Figure 6. The effect of a unit increase in flow in the clockwise direction around the 
cycle is determined by adding together the costs of the arcs whose flow increases and subtracting 
the cost of the arcs whose flow decreases. The change in objective function value per unit change in 
flow in one direction will be the negative of the change in the opposite direction. An example is 
shown in Figure 6, with the initial flow regime on the left, and the flow after cycle removal on the 
right. The supplies and demands are shown in the boxes and circles, while the numbers in 
parentheses above each arc give the cost and flow for that arc. For example, the cost C34 is 5, while 


72 


Ullhi, 



Figure 6: An initial solution with a cycle (left), and the improved flow regime after cycle removal 
(right). The numbers in parentheses above each arc are (cyjXy). 


there is initially 3 units of flow on that arc. A unit increase in flow clockwise around the cycle 
(beginning at Si) will cause a change in the objective function value of 4-64-3-8 + 5- 6 = —4, 
a net decrease. A unit increase in the counter-clockwise direction therefore yields 6, a net increase 
in the objective function. Clearly, increasing the clockwise flow is profitable, so flow is increased in 
this direction. Flow values will thus be increased for x 13 , x 2 2 and X34, while flow is decreased for 
X 23 , X32, and X14. That is, flow is increased on the arcs in the cycle which point in the profitable 
direction, and decreased on the the other arcs, until one of the decreasing arcs reaches zero flow. At 
this point, the cycle has been removed, and the value of the objective function has been decreased. 
In Figure 6, increasing the flow clockwise around the cycle by four units breaks the cycle by 
eliminating flow along X14, and reduces the value of the objective function by 16 units. The 
improved flow regime is shown on the right side of the figure. 

This te chni que is used as a local relaxation method by applying it to pairs of subproblems. Two 
subproblems which are adjacent in cost space are joined to form an M x 4 problem, which is 
inspected for cycles. If any are found, they are removed and the problem is searched again. 

Two different methods for applying this technique are investigated. The first, termed total 
relaxation, is to join adjacent pairs of M x 2 problems, remove the cycles, then repeat the process 
by joining adjacent M x 4 pairs, removing the cycles, then to join M x 8 problems, and so on, until 
the global problem for the current level is inspected and certified to be cycle- free. This approach is, 
however, extremely expensive. The second approach only employs a local relaxation, and is therfore 
true to multigrid principles. In this case, only pairs of M x 2 are checked for cycles. The gain in 
speed from using this second method is significant, while the decrease in accuracy is negligible (see 
Table 1 in the next section). 


EXPERIMENTAL RESULTS AND CONCLUSIONS 


The algorithm employed in this work is an FMG algorithm, using demand-weighting as the 
restriction method for computing costs, interpolation by local optimization, and local relaxation by 
cycle removal. The results are displayed in Table 1. While the multilevel algorithm performed well 
on problems with only two or three supply nodes, the results for the five supply node problem are 
unsatisfactory. The table clearly indicates that relaxation by cycle removal is as effective when 
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applied over a local area as when applied globally, and the computational effort required for local 
relaxation is an order of magnitude smaller. 


Problem 

Size 

Relaxation 

Method 

Run Time 

% Above 
Optimality 

2 x 1024 

Total 

1.210835 

0.02 % 

2 x 1024 

Local 

0.131437 

0.02 % 

3 x 1024 

Total 

1.15788 

8.41 % 

3 x 1024 

Local 

0.108765 

8.41 % 

5 x 1024 

Total 

1.18411 

58.4 % 

5 x 1024 

Local 

0.106392 

59.7 % 


We note also that this algorithm is not now competitive with the state of the art in network flow 
optimization methods. No numerical data are available as a careful comparison has not been made, 
however, some rough comparisons indicate that much remains to be done before a competitive 
algorithm could be obtained. 

The most significant contribution of the current research is the removal of the requirement for a 
physical interpretation of the problem, and the dependence on a relationship between distance and 
shipping costs. By mapping the problem into cost-space, a multilevel approach can be applied to a 
much broader class of problems. Of course, there is a limit to the number of supply nodes which 
this approach can handle, due to the increasing dimensionality of the problem. However, for 
problems with few supply nodes, this approach can be helpful. We predict that further work will 
yield the result that problems which have either a very small number of supply nodes, or a 
geometrical interpretation, can be solved to within an acceptable degree of optimality using a 
multilevel approach. However, problems which do not meet either of these criteria probably cannot 
be solved with currently known multilevel methods. 


Further Research 

An algorithm analogous to the full approximation scheme (FAS) should be developed. In the 
current work, we were unable to find an effective method of extracting a cornction from the 
solution on Q 2h and applying it to the approximation on while still maintaining feasibility. 
Instead, we compute the solution on fl 2/l and use interpolation to replace the solution on fl h . Since 
a direct analog to the residual in a PDE is unknown for in an optimization problem, FAS is likely 
the method of choice, however, the difficulty mentioned above must be overcome. 

Another possibility for improving this algorithm is to begin the procedure by overlaying the 
cost-space with a regular M-dimensional grid. The first step of the restriction process would then 
be to map the demand nodes from their natural irregularly spaced positions in cost-space to the 
regular grid points. Later, the final interpolation step would be to transfer from the regular grid 
back to the original demand points. This approach overcomes a shortcoming in the current 
algorithm, which aggregates demand nodes which are closest in relative distance in cost-space, 
regardless of the absolute distance between them. In using a regular grid, a demand node on Q 2h 
would reflect only the demand at nodes a distance of 2 h or less away from it. Another important 
potential advantage is that the work on each coarser level is reduced by 2 _M , instead of by one half 
as in the current research. If the regular grid approach proves worthwhile, then it could be 
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extended to a fast adaptive composite (FAC) grid approach. In a network optimization setting, this 
might be done by overlaying a fine grid on those regions of cost-space where the density of demand 
nodes is high, and a coarser grid on the areas of low density. In this way, the flow to nodes which 
are most similar to their nearest neighbors in cost-space will receive the benefit of a finer grid 
spacing, while nodes which are naturally more distinct from their neighbors will only enter the 
problem on the coarser levels. 
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